EM算法

  EM算法是一种迭代算法,是一种用于计算包含隐变量概率模型的最大似然估计方法,或极大后验概率。EM即expectation maximization,期望最大化算法。

1. 极大似然估计

  在概率模型中,若已知事件服从的分布或者其他概率模型的参数,那么我们可以通过计算得到某事件发生的概率。而在估计中,这些变成了方向过程:已知一组数据发生的结果,相当于获得了经验概率,通过这组数据假设模型服从什么分布,再通过经验概率求解模型参数。
  比如统计学校学生身高服从的概率分布,抽样1000人得到他们的身高数据,再假设身高服从正态分布,于是只需要求解$u,\sigma$即可。剩下的问题就是如何确定这两个参数,极大似然估计的思想就是求解在这两个参数下,得到的 这组样本发生的概率最大的情况便是这两个参数最有可能的取值。而抽取每个人的身高都是独立的,让每个p最大化就是让其乘积最大化,于是得到最大似然函数

  再求解其关于$\theta$的极大化问题便能得到一个对$\theta$的估计

2.EM算法

  上面是模型变量都是可观测的情况,这种情况可以直接根据观测到的样本数据,使用极大似然估计对模型参数进行估计。但若模型中含有隐变量的时候,就不能简单的使用极大似然估计进行计算。EM算法就是针对含有隐变量的概率模型参数的极大似然估计方法。
例子:
  假设三枚硬币,记为A、B、C,它们正面向上的概率分别为n、p、q。实验如下:先抛A,根据A的结果选择B或者C再次抛,将这次正面向上的结果记为1,反面记为0,其中A正面向上则选择B,反面则选择C。经过n次实验后得到n个观测值,根据这n个观测值估计n、p、q的参数值。实验过程中只能观测到最后结果,并不能观察实验过程,也就是A的结果是未知的。该模型可以表示如下

  其中y是随机变量Y的观测值,表示该次结果是1或0。z代表隐变量,表示模型中未能观测到的A的结果,$\theta=(n,p,q)$是模型参数。其中Y是不完全数据,Y+Z是完全数据。
  若是没有隐变量Z,那么可直接使用对数极大似然函数估计

  加入隐变量后,对数极大似然函数变为

  求解

  按照极大似然估计中的方法,似乎应该对$\theta,z$求导然后令其为0解方程组,然后注意此时的$L(\theta,z)$函数,log内含有求和运算,若是直接求导计算十分困难,因此退而求其次,既然要极大化$L(\theta,z)$,就先找到一个它的下界,再想办法提高这个下界,同样能达到极大化L的效果。使用$Jensen$不等式对似然函数进行放缩。

  • $Jense$不等式(琴生不等式)

  凸函数:是定义在实数域上的函数,如果对于任意的实数,都有:

   则该函数是凸函数。
  当函数是凸函数的时候,$Jense$不等式的含义是函数的期望大于等于期望的函数(凹函数相反)。图下图所示

  二维情况下可用凸函数定义来解释,当一个函数是凸函数的时候,它满足

  左边其实相当于其变量x先求期望后带入函数,右边是先求函数值再计算函数值的期望,也就是

  再回到$L(\theta,z)$中来,目的是为了将对数函数中的和项去掉,便可利用$jense$不等式的性质,将期望的函数变为函数的期望。先进行放缩

  其中最后一步用到了$jense$不等式,因为对数函数是凹函数,所以不等号反了过来,$f(E(x))\geqslant E(f(x))$,此处

  并且所添加的$Q_i(x)$满足  这是根据第三类$jense$不等式的条件设定的,不同系数的加权求和期望只要满足系数之和为1就能使用$jense$不等式。
  所以得到结论,$log\frac {P(y,z|\theta)} {Q_i(z)}$的加权平均就是$L(\theta,z)$的一个下界。这便是EM算法中E(期望)的来由。
  目前$Q(z)$还是未知的,需要根据一些条件来选择一个合适的函数,再次强调最终目的是极大化似然函数,现在我们得到了似然函数的一个下界,一个想法就是让这个下界去更好的逼近似然函数的真实值,下界是根据不等式放缩后得到的,那么当不等式取等号的时候便是最接近原函数的时候。回忆$jense$不等式$f(E(x)) \leqslant E(f(x))$,显然当$x$为常数的时候不等式成立。即

  既然要确定的是$Q_i(z)$,不妨设$\theta$为常数,事实上,EM是一种迭代算法,也就是我们是先给出$\theta$的假设值后再进行迭代计算的,设上一次为第i次迭代,$\theta ^{(i)}$为第i次的参数值,由上式得

  将c带入x

  第二步用到边缘概率,第三步条件概率。至此,我们得出了$Q(z)$的表达式,$Q(z)$是已知样本观测和模型参数下的隐变量的分布。将$Q_i(z)$视作常数去掉后,得到下面表达式

  到这儿已经完成了E步,对对数似然函数下界的期望的表示,这个期望是$logP(y,z|\theta)$关于$Q_i(z)=P(z|y,\theta^{(i)})$的期望,接下来需要做的便是极大化这个期望,也就是M步。极大化的方法就是我们常见的优化问题,使用拉格朗日乘数法即可,添加约束条件$\sum_zP(z|y,\theta ^{(i)})=1$后对$\theta,z,\lambda$求偏导令为0。求解

  迭代的停止条件常常为$||\theta^{(i+1)}-\theta^{(i)}|| \leqslant \varepsilon ,or,||L(\theta^{(i+1)})-L(\theta^{(i)}|| \leqslant \varepsilon$

3.EM算法的收敛性

  观测数据的似然函数是$L(\theta)=P(y|\theta)$,我们先证明似然函数关于每次迭代后的$\theta$是单调递增的
proof:

取对数

由前面结论

于是,对数似然函数可写作

那么

只需证右边非负,因$\theta^{(i+1)}$使原似然函数$L(\theta^{(i+1)},z)$达到极大,所以

对于第二项

  因此$logP(Y|\theta)$关于迭代的$\theta$单调递增,那么如果$P(Y|\theta)$有上界,经过$\theta$的单增迭代一定会收敛到某一值,还有一个定理是EM算法得到的$\theta$收敛值是似然函数的稳定点,人生苦短,证明就免了吧。

4.EM算法在高斯混合模型中的参数估计

  通过前面的介绍已经知道EM算法用于解决含隐变量的概率模型,EM算法可以认为是一种解决问题的思想,在不知道隐变量发生情况的条件下,我们通过计算已知观测数据,求出不同数据在隐变量的某个事件下发生的期望,再通过这个期望去最大化该期望下的似然函数,如此迭代下去便能得到一个局部极大的结果。
  高斯混合模型也是这样的一个概率模型,它由多个正态分布加权组合而成。也就是可以观测到混合模型产生的样本的分布,但是不知道这每一个样本是由其中哪个高斯模型产生的。因此,隐变量就是哪个样本选择了哪个模型。

  • 高斯混合模型

  假设观测数据$y_1,y_2,…,y_n$由高斯混合模型产生

  其中$\theta = (\alpha_1,…,\alpha_k;\theta_1,…,\theta_k)$,由前分析得隐变量是每个样本来自哪个模型,j个样本和k个模型,隐变量$\gamma_{jk}$表示如下

  于是有$j$个观测数据$jk$个位观测数据,总共$(j+1)k$个数据便组成了完全数据,完全数据的似然函数为

  其中$n_k = \sum_{j=1}^N\gamma_{jk},\sum_{k=1}^Kn_k = N$,

  显然$n_k$表示的是第k个模型出现的样本数量,因此所有模型出现的样本数量和就是N,完全数据的对数似然函数为

  确定好完全数据的对数似然函数后就可以直接套用EM算法的框架了,首先是E-step,求解对数似然函数关于隐变量的后验概率的期望

  注意其中的$n_k$是关于$\gamma_{jk}$的函数$n_k=\sum_j\gamma_{jk}$,注意上面两步的化简中,并没有直接计算期望,而是根据期望的性质对式子进行化简,将除了$\gamma$的式子看作常数,由$E(aX)=aE(x)$等期望的运算法则化简到只需要求$E(\gamma_{jk})$。
  计算$E(\gamma_{jk})$,记为$\widehat \gamma$,注意这个期望是对$P(\gamma_{jk}|y,\theta)$的期望。

  带入$L(\theta,\gamma)$得

  接下就是M步,只需最大化这个L即可

  求解方式同样是拉格朗日,约束条件$\sum_k\alpha_k = 1$